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Abstract 

We investigate the electrostatic equilibria of N discrete charges of size 1 /N on a two 
dimensional conductor (domain) . We study the distribution of the charges on symmetric 
domains including the ellipse, the hypotrochoid and various regular polygons, with an 
emphasis on understanding the distributions of the charges, as the shape of the underlying 
conductor becomes singular. We find that there are two regimes of behavior, a symmetric 
regime for smooth conductors, and a symmetry broken regime for "singular" domains. 

For smooth conductors, the locations of the charges can be determined, to within 
0{^/logN /N'^) by an integral equation due to Pommerenke [Math. Ann., 179:212-218, 
(1969)]. We present a derivation of a related (but different) integral equation, which has 
the same solutions. We also solve the equation to obtain (asymptotic) solutions which 
show universal behavior in the distribution of the charges in conductors with somewhat 
smooth cusps. 

Conductors with sharp cusps and singularities show qualitatively different behavior, 
where the symmetry of the problem is broken, and the distribution of the discrete charges 
does not respect the symmetry of the underlying domain. We investigate the symmetry 
breaking both theoretically, and numerically, and find good agreement between our theory 
and the numerics. We also find that the universality in the distribution of the charges 
near the cusps persists in the symmetry broken regime, although this distribution is very 
different from the one given by the integral equation. 

keywords: Fekete points; singular domains; symmetry breaking; electrostatic equilibrium. 

1 Introduction 

The placement of charges on a conductor is a classical problem in both physics and mathe 
matics. In physics the history of the problem goes back to the work of J.J. Thomson |351 [7| l6] 

* mkb@uchicago.edu 

^ isabelle@control.uchicago.edu 

*leop@uchicago.edu 

^ nicewicz@uchicago.edu 

' scvenkat ©control . uchicago . edu 



1 



The placement of charges on a hne was investigated by Stieltjes |33j in the late 19th century. 
Fekete recognized the connection between the placement of charges on a 2 dimensional 
domain, and questions involving polynomial functions of a complex variable. This problem 
was also investigated by Frostman |14j, for its connections to Potential theory. In general, the 
problem is to place charges of strength \/N on a surface so that the energy is minimized. 
In this placement the force component parallel to the surface for each charge vanishes. 

This problem has quite a different nature depending upon the dimension of the system. 
In three dimensional space, one uses an inverse square force law. The placement upon any 
surface, even one as simple as the surface of a sphere is quite complex. The discrete charges 
form an approximation of a lattice with many defects |lH [7|l6ll8]. The solution might properly 
be described as being chaotic. 

In contrast, in two dimensions, one finds smoother behavior. The problem is one of 
placement of lines of charge on the curve which bounds a two- dimensional conductor. Forces 
between the charges are given by an inverse first power of the distance. 

The analogous continuum problem has a very good general formulation via the Riemann 
mapping theorem. To calculate the distribution of continuum charges on the exterior of 
a simply connected region in the plane, construct the unique function F{w) that takes the 
exterior of the unit circle into the exterior of the region, and has the property that F{w) w 
as w ^ oo. Then the density of charges at the point z = F{w) = F{e^^) on the surface of 
the conductor is given by 

p{e) = \wF'{w)\-^ (1) 

This result provides the basis for all subsequent work on discrete charges. 

More recent work |221 1231 124j has looked at the placement of discrete charges in terms of 
the placement of points on the circle at 

Wj = e^"^^^ (2) 

The 0j's are real. According to a set of important recent theorems |^, for sufficiently smooth 
curves C , we can determine the charge placement for all large enough values of N in terms 
of a smooth periodic function, tpiOy, by writing 



9 + + ^ 



(3) 

With appropriate choice of the expression |21 gives the placement of charges with error of 
order ylog(iV)/A^^ [25] for all extrema of the energy 

Here we represent the usual logarithmic potential for line-charges by taking the real part of 
the complex logarithm function. 



1.1 New empirical results 

We calculate the actual position of charges on two different kinds of shapes: ellipses and 
hypotrochoids. These two kinds of shapes are defined as the image of the unit circle under 
maps of the form 

Fiw) = w-^ (5) 
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Figure 1: Placement of 16 charges on two different ellipses, c = 0.2 and c = 0.7 respectively. 
Notice the breaking of the parity symmetry in parts b. The ellipse degenerates to a line when 
c = 1. 



Here r = 1 for ellipses and r = 2 for hypotrochoids. The parameter c controls how different 
these figures are from the circle. 

Results of these mappings are shown in figures [TO These figures also include the positions 
of the charges placed upon these conductors. As one might expect from equation^ the charge 
densities are highest in the high curvature areas of the bounding curves. However, the figures 
also contain some surprises. In both cases, there is a symmetry breaking in the more pointed 
figures. In this breaking, the left-right symmetry is broken at each pointy region, and the 
overall parity symmetry of the figure is lost. We shall explore this symmetry-breaking in 
more detail below. 

To understand these charge placements, look at the energies they generate. When the 
sizes of the figures are properly adjusted the continuum correlation energy is zero. The next 
term in the energy is the self-energy: 

This term is independent of charge placement. The remaining energy is a correlation energy, 
and has a power-law behavior in N for large A^: 

~ (7) 

Here the exponent s is known to be two for smooth curves. In figureOlwe plot the correlation 
energy versus for three different hypotrochoids. In the first case, we can clearly that s = 2. 
For the singular figure, corresponding to figureEfc, the log-log plot exhibits a slope of 1.5. In 
the intermediate case, corresponding to the blunted points of figure ^p, there is a crossover 
from s = 1.5 for smaller A^ to s = 2 for larger A^. 

^This functional is conventionally denoted by "1>. We however reserve $ for the electrostatic potential 
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Figure 2: Placement of 18 charges on three different hypotrochoids, corresponding to c = 0.2, 
c = 0.45 and c = 0.5 respectively. Notice the breaking of the parity symmetry in part b and 
c. The hypotrochoid is smooth for < c < 1/2 but gains three pointed cusps at c = 1/2. 
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Figure 3: Plot of correlation energy versus N. The three curves shown correspond to the 
three hypotrochoids shown in figure |21 
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1.2 Connection to Modern work in dynamical systems theory 

Note the close connection between equation |31 and the expression giving the positions of 
points in long cycles and quasi-periodic orbits of KAM [5j theory. KAM theory provides 
an expansion about continuum motion, in which the positioning of points in a cycle closely 
follows that of the corresponding continuum orbit. In this Coulomb system too, a near 
continuum result gives a charge density which is only a little different from the one in the 
continuum case. As in KAM theory, equation |31 describes a near-continuum result by giving 
a function, here ip, which describes the deviation from and permits a smooth passage to the 
continuum limit in which goes to infinity. 

In KAM theory, the cycles closely follow a continuum orbit. So the continuum theory 
predicts the curve on which the points will occur; but it does not give the positioning of the 
points in the orbit. A more subtle calculation is required to find that. Here the analogous 
fact is that the initial angle 9q of equation |21 is not determined by any continuum theory. In 
fact, one does not know a first principles method of finding ^o- 

Dynamical systems theory is much concerned with bifurcations, i.e., situations in which 
a qualitative change in behavior occurs. In the Coulomb problem there is a bifurcation from 
the single symmetric to multiple solutions and symmetry breaking behavior. In dynamical 
systems theory when a parameter is varied, bifurcations can occur for longer and longer orbits 
producing an interesting kind of critical behavior jl7 ( I2()| 1^. Here the analogous situation 
arises when the bounding curve develops cusps as the parameter c is varied, and close to 
this cusp formation, the Coulomb system exhibits a symmetry breaking bifurcation. This 
behavior produces interesting anomalies including the cross-over in the A'^-dependence of the 
energy. 

1.3 Paper Outline 

Section |2 of this paper deals with the equations for the placement of charges on the circle 
and on the circle mapping onto "real space". We find that, for smooth curves, the uniform 
placement on the circle gives a correlation energy which goes exponentially to zero for large- 
N. We derive an integral equation which describes the "continuum" behavior of the deviation 
of the charges from the uniform placement on the circle. This is related to, but different from, 
the equation formulated by Pommerenke [29l. Our derivation gives a physical interpretation 
of the equation as being related to the distortions of fields produced by the discrete charges' 
self-energy effects. 

Section |S1 follows the analysis of the integral equation in Fourier space. Here the first 
quantity of interest is the structure function Sk, which is the Fourier transform of the charge 
distribution. The equation is solved in the large A- limit for both the ellipse (an old result ^2i^ ) 
and the hypotrochoid-a new result. The solution, involving an elliptic function, is formally 
very similar in the two cases. We also compare the solutions we obtain with numerical 
observations, and we find that the Integral equation is not applicable for "singular" shapes, 
and in the symmetry broken regime. 

In Sec. HJ we present new results for the singular limit of the Coulomb problem. We 
investigate the energy and the charge placement in the symmetry broken regime, i.e., for 
ellipses with c ~ 1, hypotrochoids with c ~ ^. We also study "singular" curves with cusps 
and corners, viz. the line segment, the hypotrochoid with c = ^ and various regular polygons. 
We present scaling arguments for the dependence of the energy on N in the singular limit, 
and use this to obtain scaling laws for the symmetry-breaking transition. We also compare 
these results with numerical simulation. 



6 



2 Real Space Analysis 



^(z) = -— \n{z - Zj) 



2.1 On the circle 

We describe our electrostatics problem by using complex variables. Thus a charge of size 
1/N placed at the point {xj,yj) will generate the complex potential at the point {x,y) via 
the formula 

1 

N 

where In is the complex logarithm and z is x + iy. If there are N such charges 

1 ^ 

Hz) = -^Y.^n{z-z,) (8) 

i=i 

When the charges are uniformly distributed on a circle of radius one, then the potential is 

1 ^ 

«>(z) = -— J^ln(z-e2^^j/^) 

= -lln(z^-l) (9) 

Note the elegance of the result of equation (jUJ. In the limit of large N, inside the unit 
circle z'^ is very small, and the potential is well approximated by $ = {—iri)/N. Outside, 
we find ^ = —Inz which is the continuum result. To get the potential seen by the charge at 
z = 1, we have to subtract away the charge potential it produces and use 

^°W = ^i-;^ (10) 

Equation (fTU]) describes the potential upon our charge as the potential produced by the 
charges on the circle less the potential generated by our charge itself. In the neighborhood 
of 2; = 1, we can evaluate this potential as 

The total energy of this set of charges is obtain as one half this potential (at z = 1) times 
the charge summed over all charges. The result is 

log 

Ej^ = -'^^ = E, (11) 

We previously described this terms as a self-energy. More precisely the energy of our charges 
is the energy generated by a continuous circle of charge, namely zero, minus the energy needed 
to put together the "missing" portions of the circle. Those portions are the lines of charge 
counted as not present in expressions like equation (|1U() . These are lines of charge of length 
2Tr/N not generating potential in the gap around our charge itself. The energy of one such 
gap is (log A^)/2A^^. Subtracting away N such energies gives expression ((TT|) . 

2.2 Uniform Spacing on the circle 

This polynomial method is elegant but it is difficult to generalize to other curves beyond 
the circle. A method less elegant but more susceptible to generalization is to evaluate the 
potential via the contour integral: 
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Figure 4: Contour for calculating the potential generated by — 1 charges on the circle. 



Here the contour, C, encircles all the poles of l/(y^ — 1) except the pole at y = 1 (See Figure 
0}. By summing the contributions from the various poles we find the result of equation Q 
once more. 

An alternative mode of evaluation is to deform the contour to include the poles at zero 
and one and the branch line which runs from z to infinity. The result is then 

og.+M^-ir^^ (13) 

N N J s s — \ 

which works out to the result Q once more. 

The next step is to generalize the calculation to include the possibility of a curve other 
than the unit circle. Let F{w) be the function which conformally maps the exterior of the 
unit circle into the exterior of a curve C. We shall place our charges on the curve C, but 
give them the placement which is appropriate for the mapping of charges uniformly spaced 
on the unit circle. Now consider once more the effect of — 1 charges on C. Instead of the 
integral (|12j) . the potential on a given charge is set by the more general object 

^n{z, v) = - — ^ log{z - F{w)) (14) 

Here v is on the unit circle and the contour, C, encircles the points 

wj =ve'^'''^/^ for j = 1,2,... ,iV- 1. (15) 

in the positive sense (See figure El) The real part of $Ar(z,t;) is the potential at the space 
point z, produced by the specified set of charges. Note that this potential does not take into 
account the redistribution of charges required to produce an equilibrium distribution. 

The calculation follows much the same lines as the previous one, but the details of the 
calculation are a bit different because the integrand is different. As we move the contour we 
get four different contributions to the result (|14)) (See Figure IHI) The first contribution comes 
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Figure 5: The contour for equation (|14() . The main difference from the previous figure is 
that there are some branch hues within the unit circle. This figure is drawn with v = 1. 

from the encirchng of all the branch lines of the logarithm. The number of the branch lines 
that connect the preimages of z to depends on the map. The next contribution comes from 
encircling the pole at w = v. This gives a contribution from the "missing" charge, which is 

^log[z-F{v)] 

Next comes a contribution from circling the branch line which goes to the point at infinity. 
Because that encirclement contribution has a logarithmic divergence, we do not carry it out 
all the way to infinity, but instead cut it off at some large but finite length L ^ 1. Then 
we also need a contribution from a circle at infinity. All this is illustrated in figure El The 
branch line contribution is 

dw 1 f^^ du u 



F-^{z) ^ — v'^ N J^p-if^^-^^N u u — 

1 , L^-t;^ 



1 



\ogL--\og[[F-\z)] 



N 



Finally, the large circle gives just 

-logL 

since for large w, F{w) is approximately w. When we add up all these terms we see that 

^^{z,v) = ^\og{z-F{v))-^ log(n^-t;^) (16) 

MeF-i{z} 

Equation 1)161) is a most interesting result. For z inside the curve C, the potential is 

- log(-f ) + log[z - F{v)] 
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Figure 6: The contour for equation (|14l) after deformation. 



The first term is an imaginary constant describing how we placed the branch hne. The second 
is a contribution of the "missing" charge. For z outside the the curve C the potential is 

-log {F-\z)\^^\o^{z-F{v)\ 

The first term defines the potential in the limit in which N goes to infinity. This is the familiar 
but deep statement that the inverse of the mapping is proportional to the exponential of the 
complex potential. The second term is once more the missing charge. 

Finally, we can look at the case in which the potential is evaluated at the unrelaxed 
position of the iVth charge, i.e. z = F{v). We then have: 



^n{F{v),v) = lim ^n{z,v) = lim ^ log 



z - F{v) 



n 



mGF-1{z} 



Clearly, it is f G ^{F{v)}, so all but one part of the limit can be trivially evaluated. The 
one piece in which u = v can be done by I'Hospital's rule. It then follows that: 



^NiF{v),v) = ^log 



F'jv) 



- T 



log(n^-t;^^) 



uG-F-i{F{v)},u^-!; 

Now we can write down the total energy of the system. Since 

N 

2N 

we have 



1 ^ 

EN = ^Y.^N{F{vj),V,) 



(17) 
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In general, it is difficult to find all such u' s and perform the sum. However, in the case 
of the ellipse given by F{z) = 2 + |, this expression nicely simplifies. 

2.3 Ellipse energy - uniform spacing 

To begin, let us cvahiatc the first sum in the above expression for the uniform spacing of 
charges, i.e. put Vj = Q, where Q's are the iV*'* roots of unity. 



We have. 



N N 

^Emc,^'(0)) = ^ log 11(0- 
j=i j=i 



It follows that 



N N 

C N / AT "TT / ^ C 



11(0-^) = (-1)^11(1- 

j=l ^1 j=l '^J 



We notice that since the roots of unity form the cyclic group of order N, and taking N to 
be odd, 2 always has a multiplicative inverse in Zjv, and the group is invariant under squaring 
the roots of unity. 

Therefore, for odd N, 

N N N N 

i-^f 11(1 - ^) = - 11(1 - ^) = - 11(1 - cC^^-,) = -c^ n(^ - CiV-.) = 
j=i j=i j=i j=i 

= -(l-c^) 

It follows ^ 

2^ E log(0^'(0)) = 2^ log(l - c^) 

If N is even, then 2 * Zjv = Zjv/2 so that 

N N 



(-1)"" 11(1 -72)= n (1 - ;^)^ = (1 - c^/^)^ 

and therefore 



j=l i=2,4,6... 



1 ^ 1 

^j:^og{CjF%)) = ^ log(l - c^/^) 

i=i 

Please note that the evaluation of this sum can be readily generalized to the family of 
curves given by F{z) = z + ^(p-i) where p is prime. Then zF'{z) = z ~ ^^(pll) , and of course, 
as before, ii N ^ mod p, raising the A''*'* roots of unity to the p*'^ power simply reshuffles 
them so that 
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1 ^ 1 

2]^El°g(0^'(0)) = ^log(l - (b- l)c)^^ 



In the case that = mod p, 

11(1 -^)= n a-^f-d-ap-i).)*)' 

j=l S j=p,2p,3p... 

and we conclude that 

1 ^ 

^ j;iog(Or(0)) = ^log(l - ((p- l)c)f ). 
i=i 

Now we just need to evaluate the second sum for the ellipse and we will have the exact 
energy expression for a uniform arrangement of charges. To this end we need to obtain the 
u's which satisfy the equation -F(Cj) = Fiu) or — uF{(j) + c = 0. It is clear that u = (j 
satisfies the equation, and from this it follows that, since the product of the two roots of the 
above quadratic is c, the other u must be -r-. 

We now have: 

1-,y: e iog(""-"'')=4iog(i-c'') 



2iV2 ^ ^ ' 2N 

i=l uGF-^{F{v)},Uy^v 



since =1. 



Finally, for odd N, 

logiV 1 



Ej^ = _ ^ iog(i _ c^) + -i- iog(i - c^) 

while for even N, we have: 

The remaining material in this chapter should be considered to be more heuristic in 
content. 



2.4 Heuristic Discussion of Energies 

We have just completed a careful analysis of the self-energy in the case in which the charges 
are uniformly distributed on the unit circle. Now we extend this discussion to include the 
effects of a non-constant charge density on the circle. In fact, if we define the charge density 
to be the smooth function of j defined so that 

ej-eK= f ^ (18) 

If pj is independent of j we know that equation 1)11(1 gives an essentially exact answer for the 
large-iV case. 

Our result for contributions of the self-energy to the potential for a single particle is given 
by equation (|17() . which can be written in terms of the spacing between particles on the curve 
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C which is given by a spacing factor on the circle times the conversion factor from circle it 
curve \/F' . Thus the spacing is 

2NTip ^ ' 

If we take our original expression for the self-energy contribution to the potential, and then 
say that this contribution to the real potential depends only on the actual spacing we find a 
potential contribution 

1 - r|F'(^;,0^ 



N 



log 



(20) 



Just as before, we calculate a contribution to the total energy by taking half the potential 
multiplied by the charge density and adding. We find 



1 ^ 
2iV^^° 



\F'{v,) 



2NTTpj 



(21) 



This energy is the continuum value (zero) plus a contribution from the replacement of 
the line of charge distributions immediately surrounding a given charge (which is included 
in the continuum energy) by the self energy of the point charge which replaces it, and is not 
included in the energetic calculation. This is the charge in question 1/N times the separation 
of the charges on the line, logl/A^. One can also write this energy as an integral over the 
circle of the form 



1 

2N 



2tt 



de p{9) log 



\F'{v)\ 

2NTTp{e) 



(22) 



Here v is e*^ and p{9j) is the same as pj. This form of the energy is convenient for variational 
calculations. 

Note the factor of one half in front of equation (|22() . It will be important in what follows. 
The same discussion can be easily extended to the calculation of the interaction energy. 
Because of the logarithmic potential the result is 



Er. 



2tt 



de p{d) 



2-K 



dp p{p) 3f?log[F(e^^) -F(e^^)] 



(23) 



2.5 Variational Calculation 

To construct a variational principle, we take the two energies of equations (|22|) and (|23|) add 
to them a Lagrange multiplier term of the form 

de p{e)<Po 







designed to ensure that the total change is constrained to be unity. One then sets the resulting 
variation with respect to p to zero and finds 



j^J dp p{p) Klog[F(e^^) - F(e*^)] + ^ log 



2N-Kep{e) ^ ' 



The e in the last term appears as a result of varying the density inside the logarithm. Note 
once more the factor 1/2 in front of the last term. This equation is as far as we can tell, new. 

Equation (|24j) looks like the statement that the electrical potential on a conductor is 
constant. However, that is not quite the right interpretation. To the requisite order (1/A^) 
the right hand side of equation ([2i|l is not the electrical potential. The electric potential is 
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the somewhat shnilar expression in which the symbols e and 1/2 just mentioned are replaced 
by unity, i.e., 



2tt 



1 



dp p{p) mog[F{e'") - F(e*^)] + - log 



|F'(e* 



2N'Kp{e)_ 



This actual potential, in contrast to varies in space. 

In equation (|24j) the potential </>o can be evaluated by integrating the equation over all 
values of 6. The very same arguments which we used in section 1221 to perform the integration 
over charges uniformly distributed on the circle can be used to simplify the expression for (pQ . 
We find 



27r 



de — log 

2N ^ 



\F'{e 



2N7Tep{e) 



We are looking for the lowest order terms which are of order 1/A^. Hence we can set p{6) 
to its lowest order value l/(27r). We can also evaluate the integral over logF' by taking the 
contour to infinity. In the end we find 



27r(/)o 



27r 



1 



Ne 



so that, to order 1/A^ equation (|24() becomes 







2tt 



dp p{p) mog[F{e 



-F{e^>^)] + —\og\F'{e'') 



(25) 



which is our desired final result. 

This equation is related to the integral equation for the charge distribution in Ref. |29j . 
While our equation is a statement of the constancy of an appropriate potential on the surface 
of the conductor, the equation in Ref. j29j is a statement that the tangential component of 
the electric field is zero at each charge location, in an equilibrium configuration. 



3 Fourier space analysis 

We will analyze the integral equation (|25|) using Fourier analysis techniques. We define the 
Fourier coefficients of the continuum density p by 

Pk= p{p)e'^'^dp 
Jo 

Following Pommerenke jJH], we also define 

N 

where Wi is the location of the ith charge in the minimum energy configuration. The 
encode information about the locations of the charges. For k <^ N, the Sk reflect the smooth, 
large scale behavior of the distribution of the charges. Thus, they should correspond to the 

Pk by 

/•27r 

Sk^Npk = N / p{p)e'^^'dp. 
Jo 

For larger k, in particular for A; ~ A^, the discreteness of the charges become important, and 
Sk will in general be very different from pk- 
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We now rewrite the integral equation ()25() in terms of these definitions. From H25|l we 
obtain 



27r 



\og[e'^ - e^^] + log 



Substituting p{^) = ^ p^e we obtain 



Pk , sr^ - Ck 

^+l^^MPl = ^, 

I 

where, following Pommerenke |2^, we have defined the matrix aki and the vector Cj by 



k=l 1=1 



and Cj = - Y^k=i (^k,j-k, so that 



logF'H = -J2 



CkW 



k=l 



Solving this system of linear equations will yield the continuum density p. Note that the 
above system is only meaningful for A; = 1, 2, 3, . . .. We know that the total charge is 1, and 
this fixes po = 1. We can also obtain information on the locations of the charges form the 
correspondence between pk and Sk- For large A^, using the correspondence Sk ~ Npk, we get 



k 



akisi 



Ck 

2 '■ 



(26) 



which is the same as the equation in Lemma 2 of Ref. [2H1- We will henceforth refer to this 
as the Pommerenke equation. Note that our analysis yields the (unexpected) conclusion that 
the equations determining the Sk have no dependence on A^, the total number of charges. 

We can calculate Uki and Ck for the ellipse and the hypotrochoid using the explicit forms 
of F. In order to keep the formulae simple, and avoid the appearance of complex roots of 
unity, we will use the forms 

FiO=C + -^, 

Therefore, the c in this section is the same as the c in the previous section, if we set ^ = e™^'^w 
for the ellipse, and ^ = e^^w for the hypotrochoid. We thus need to shift the results with our 
choice of mapping function by 7r/2 and tt respectively to compare with the mapping functions 
in Sec. [21 

A calculation for the ellipse yields 



Ckl 



and 



Ck 



2 k/2 
I, <^ 



k odd 
k even 



(27) 



(28) 



A similar calculation for the hypotrochoids shows that the energy can be put in a similar 
form with 

k + l^O mod 3 

k <p or I <p,p= {k + l)/3 (29) 

f {k-pW-py- otherwise, p = (/c + /)/3 
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Figure 7: as a function of k for the ellipse, for c=0.2, 0.6, 0.8 and 1 respectively; the circles 
are the numerics obtained with N = 40, the stars with = 100, while the dashed line is 
expression 1311 As expected the agreement is good for small values of c, when the symmetry 
is preserved. The deviations are therefore stronger for = 40 than for N = 100. 



and 



/ mod 3 

|(2c)'^/^ otherwise. 



(30) 



3.1 Solving the Pommerenke equation 

For the ellipse, the matrix a^i is diagonal. Consequently, Eq. (|26() is easily solved to yield: 







Sk 



1+c* 



k odd, 
k even. 



(31) 



which is the solution presented in 28 . Note that for this solution, Sk is real for all k. This is 
equivalent to the statement that there is no symmetry breaking. Thus this approach cannot 
capture the symmetry breaking for the ellipse, that we see in the numerics. 

This result is compared to numerical calculations in Fig. [7| which is also done with the 
map F{w) = w + c/w. The circles are the results obtained for N = 40, while the stars 
correspond to = 100. As expected, the agreement is good for small enough values of c, 
for which the symmetry is unbroken. Moreover the agreement is better for = 100 than for 
A^ = 40. This is in agreement with the expectation that deviation of the solution 1)311) from 
the true Sk is 0{k'^/N) 

We will now consider the hypotrochoid. We first enumerate the consequences of the 
three-fold symmetry as well as the reflection symmetry of the hypotrochoid. We get, 

1. Sfc = for k ^ mod 3 is a reflection of the symmetry of the underlying hypotrochoid 
under rotations by 27r/3. This symmetry is preserved in our numerically obtained 
solutions if the number of charges A^ is a multiple of 3, and is broken otherwise. 
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2. Sk is real as a consequence of the symmetry of the hypotrochoid under reflections about 
the real axis. This symmetry is broken in our numerical solutions for any N, provided 
that c is sufficiently close to 1/2. 

If there is a unique minimizing configuration for the electrostatic energy of the hypotrochoid, 
the symmetries of the problem imply that Sfc = if A: 7^ mod 3 and Sk is real for all k. 
These features are also reflected in the (approximate) obtained by solving the Pommerenke 
equation. Since the Pommerenke equation is independent of A^, we cannot investigate the 
dependent features, like symmetry breaking for the hypotrochoid, using this approach. 

We will now solve the Pommerenke equation for the hypotrochoid. If c < 1/2, and 
the aki decay exponentially as k and I get large. Therefore, it is reasonable to try the series 
solution 



kck sr-fi , 
Sfc = ^ - l^ikaki)— + 



^^{kaki){laij)-^ + 
3 I 



This series converges (not just for the hypotrochoid, but for any admissible ati)^ as proved 
by Pommerenke j28j . 

Let s^j^ denote the i-th term in the series, so that 



kck 
2 



y^^kau 



'Z ' 



2,3, 



Eqs. (HH) and together give 



0, ky^O mod 3 



:(2c) 



k 

E 



kcP 

p 



{k - p)\{2p - k)\ 



1,2,3, 



(32) 



Since 3p - 
Setting k ■ 



A; 7^ mod 3 unless A; = mod 3, it follows that = if A; 7^ mod 3 for all i. 
= 3m, p = m + n and replacing the sum over p by a sum over n = p — m, we get 



2m 

E 



[ail] 



f 



2m 

E 

_rm+ii 



(m + n)l 



3m 



m + n J (2m — n)!(2n — m)l 



m 



n + m 



m + n 
2m — n 



,3 \ ™+" 




^3n 



1,2,3,... 



We can evaluate these sums using the saddle point method. Using Lemma ^ in Ap- 
pendix and = fr-Q™, we obtain the asymptotic expressions 



3m 



1,2,3,... 



where 



f3i 

Pi+l 



^0 



64A 



1,2,3,... 
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and all corrections to this asymptotic result are 0{l/m). 

At the critical parameter value c = ^, ro = 1. In this case, /3i = 1 and 



•i+l 



64A 



implies that Pi = 1 for all i. Thus, the series solution breaks down in the critical case, as the 
solutions sj^ do not decay with i. 

For parameter values close to the critical value, c = ^e"*^, we have 



= 2c = e"' PS 1 - e 

If e is sufficiently small, then we have Pi is close to one for the first few values of i. We can 
study the dependence of Pi on the parameter e by setting Pi = e~*»^+*^('^ ) = l — t^e + O(e^). 
Linearizing the relation between /Jj+i and Pi close to the fixed point /? = l,ro = 1, we get 

Pi+i = l-iU + 2)e + 0(e2) + . . . 

This gives 

ti+l = ti + 2. 

Combining this with = 1, we get 

ti = {2i - 1) 

It is clear that this analysis is valid long as <C 1, i.e., for i < e^^. Beyond this value, we 
can show that the Pi decay rapidly. Using this expression for the Pi, we can re-sum the series 
solution for to obtain 

The complete solution is 

r kj^O mod 3, 

^ I iTTi^(l + Oik-')) ^It0^ = mod 3. 

Note the similarity with the solution in Eq. ()31() . The agreement of result (|33|) with numerics 
is tested in Fig.|HJ We use the map F{w) = w — c/uP', but the quantities plotted are (— l)'^Sfc, 
which is the same thing as using the map F{w) = w + c/w^ . For small c and for k N, the 
agreement is very good. However, the agreement is rather poor if c ~ 1/2, indicating that 
the integral equation is no longer valid as we approach the singular limit c 1/2. 



3.2 The minimum energy configuration of the charges 

Our goal is to find the locations of the N charges on the ellipse and hypotrochoid that 
minimize the electrostatic energy E]\[. The results from the previous section tell us that 
the locations Wi on the unit circle that give the appropriate locations on the ellipse and 
hypotrochoid by Zi = F{wi) are such that 

k even, on ellipse. 
+ A; = mod 3, on hypotrochoid. (34) 

otherwise 
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Figure 8: Sk as a function of k for the hypotrochoid, for c=0.1, 0.3, 0.4 and 0.5 respectively. 
As expected the agreement is good for small values of c, when the symmetry is preserved. 



for 1 < k N. The problem of finding the locations Wi reduces to the problem of solving 
the above set of equations. 

As Pommerenke shows in Ref. the locations of the charges are given in terms of the 
solution Sk to the Pommerenke equation by 



9 + + ^ 



(35) 



where 



oo 



k=l 



ike _ ^ ike 
k 



and ^0 is an initial angle that is not determined by this analysis. 
For the ellipse and the hypotrochoid, we get the series 



oo ^ 
k=l 



+ 



^2k 



sm(mk6), 



where m = 2, 3 for the ellipse and hypotrochoid respectively, and q = r™. This series can be 
summed using the formula 



oo „ 

aiaiu) = h 2 > — 

n=l 



mru 
TT- sm — — . 

2n 



(See p. 277, Ryshik and Gradstein to yield 

mO 1 



-am 
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K = K(k) is the complete Jacobi elliptic function of the first kind PH], and the modulus k 
is determined by requiring that q be the nome pQ, 



exp 



TlK'ik) 

' K{k) 



where K'{k) is the complementary complete elliptic function of the first kind, K'(k) = 

We are particularly interested in situations close to the critical value of c. In this case, 
q = PS 1 — e is close to 1. From this, we infer that K' /K ~ 0, i.e k is close to 1. We will 
now work out the asymptotics close to k = l,q = 1. 

If A; = 1 — (5, A;^ « 1 - 2(5 and k' = \J\ — k'^ ~ \/25, and we have 



K{k) 



Tv/2 



dx 



\/l — k'^ sin^ X 

dx 



cos'^(x) 



= -—iog{6)+e{s) 

where we have extracted the dominant singularity of the integral in the second line, and 
are functions that remain bounded as 5 — > 0. We also have 



w/2 



dx 



'0 Vl - 2(5 sin^ X 

Combining this with the previous result, and the definition of the nome q, we see that 

^2 



q « exp 



vr 



LV21og(,5) 



Rearranging this gives the result 

^2 



V2 



log(g) log((5) ^0 as (5 ^ 0. 



Using log{q) = — e, we see that 5 ~ exp[— vr /(^/2e)], so that K log (5/^/2 ~ vr / (2e). 

Fig. ini is a comparison between the results of this calculation, and the numerically ob- 
served displacements of the charges form the positions that would correspond to a uniform 
distribution on the circle. The normalized deviation function 



ip{9) = — - -am 



vr 



(36) 



corresponds to the following two changes relative to the uniform distribution on the unit 
circle - 

1. The "charge" density in the smooth regions is increased by an amount bpj p = m/{4N). 

2. The charge moves away from the vicinity of the cusps, and this depletion region has a 
width Ti/{mK) ~ 2e/(mvr) around the the angles corresponding to the cusps. 

Note that Figs.|ni(b), (c) and (d) are obtained with the mapping F{w) = w — cjw^ and the 
solid lines are therefore given by the "shifted" function ^(0 + vr). 
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Figure 9: A comparison between the theoretical prediction and the actual locations of the 
charges in the minimum energy configuration. The normalized deviation from the nominal 
position 9j = 27r j/N is plotted, (a) shows the deviations for an ellipse with c = 0.2. (b),(c) 
and (d) depict the charge placement on a hypotrochoid with c = 0.2, 0.3 and 0.4 respectively. 
The solid lines are the theoretical prediction. 
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4 The singular limit 



In this section, we investigate various aspects of the Coulomb problem on domains that are 
either singular, i.e., have cusps or corners, or are close to being singular. We begin by first 
describing some of the numerical observations, which motivate the theoretical analysis that 
follows. 

4.1 Ellipse, hypotrochoid and related curves. 

4.1.1 Symmetry breaking 

As previously mentioned, an interesting feature of the charge placement on the ellipse and 
on the hypotrochoid is that the minimum energy configurations break the symmetry of the 
problem in certain ranges of c and N. For a given number of charges A^, if the value of the 
parameter c is larger than a given threshold c* which depends on N, the mirror symmetry 
is broken around each pointy region. The larger A^, the closer the threshold is to the critical 
value, i.e, c*{N) — > 1 for the ellipse and c*{N) — > ^ for the hypotrochoid as A^ ^ oo. 

This transition can be described by a phase diagram where c* is plotted versus A^: the 
case of the ellipse is presented in Fig. ^Jand the case of the hypotrochoid in Fig. 1111 

On both phase diagrams, two different cases should be distinguished, depending on the 
number of charges A^. Let us first consider the case of the ellipse. If A^ is even, in the unbroken 
symmetry phase, the charge configuration will be symmetric versus both the vertical and the 
horizontal axis. Both symmetries break when c > c* (see Figure^. On the contrary, for 
odd A^, there is no symmetry across the long axis, while there is always symmetry across the 
short axis (see Figure [T2|l . In this case, there is no transition of any sort. The value c* at 
which the symmetry breaks for even N forms a curve (see Figure ll()|) . 

Similarly, in the case of the hypotrochoid, the case A^ = mod 3 has to be distinguished 
from the case N ^ mod 3. In the first case, a symmetric charge configuration respects the 
three mirror symmetries, one at each pointy region (see Figure 12)). In the case A^ 7^ mod 
3, only 2 out of the 3 branches of the hypotrochoid are identical, with one charge more or 
less than the third branch. Therefore, only one mirror symmetry is respected in this case 
(see Figure [T3|) . The boundary between the broken and the unbroken symmetry phases is 
different for A^ = mod 3 and for A^ 7^ mod 3 (see Figure [TT|). 

4.1.2 Energy dependence on A^ 

• Line segment: As pointed out previously, a useful quantity to characterize the placement 
of charges on a given curve is the energy as a function of the number of charges. The 
case of the line segment r2 = [—2,2] has been treated analytically in literature |27j . 
The energy dependence on A^ is given by 

_ InN ln2 InA^ A 

where ^ is a constant. The first term is the self-energy. The first correction, or corre- 
lation energy, goes as A^~* with s = 1 in this case. 

Other curves presenting this kind of A^-dependence for the correlation energy are the 
star shape curves defined by the conformal mapping F{w) = z{l + cu;~™)^/™^27 , with 
c = 1, where w = e*^ is on the unit circle and m is the number of star branches. Figure 
1141 shows the case m = 3. The dependence of the correlation energy on A^ is shown in 
Fig. Notice that the charges are grouped near the points and that the density is 
lower at the center. 
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Figure 12: Placement of 15 charges on two different ellipses, c = 0.1 and c = 0.7 respectively. 
Notice the breaking of the parity symmetry in part b. 



c=0.3 



c=0.4 





Figure 13: Placement of 19 charges on two different hypotrochoids, c = 
respectively. Notice the breaking of the parity symmetry in part b. 
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0.3 and c = 0.4 



Figure 14: Placement of 15 charges on a star shape curve with 3 branches for c = 1. 




Figure 15: Dependence of the correlation energy on N for the star shape curve with 3 
branches. The slope shown for the data is consistent with the value -1.0. 
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• Ellipse: The line segment r2 is the singular limit, c = 1, of the ellipse defined by the 
conformal map F{w) = w — For a smooth curve, c < 1, the correlation energy goes 
as with s = 2 (see Fig. I16() . See the work of Pommeranke where it is shown 
that s=2 for any smooth curve. In the case c close to but smaller than 1, we observe 
a crossover as a function of the number of charges A^: for smaller values of A^, s = 1, 
while for larger values of A^, s = 2, as shown in Fig. 1161 




-7 -I , 1 , , 

12 3 4 

Note:E=E + ln(N)fN log(N} 



Figure 16: Dependence of the correlation energy on N for the for the ellipse, with c = 0.975 

The same dependences are observed for the star shape curves for < c < 1 (see Fig. 
[T7|) . The different situations are depicted in Fig. [THl 




Figure 17: Placement of 30 charges on a star shape curve with 3 branches, with c = 0.5. 

• Hypotrochoid: In its singular limit, c = 0.5, the hypotrochoid has a correlation energy 
of the form A^^'* with s = 1.5 according to the theory described below. The data 
supports this result (see Fig. whereas s = 2 for smooth cases c < 0.5. As for the 
previous curves, a crossover is observed from s = 1.5 for smaller A^ to s = 2 for larger 
A^ when c is close to but smaller than 0.5 (see Fig. inj. 
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0.01 




Figure 18: Dependence of the correlation energy on A'^ for the star shape curve with 3 
branches, for c = 0.5 and c = 0.95. The slopes of the line is consistent with the theoretical 
values of -1.0 and -2.0. 
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4.2 Regular polygons 



4.2.1 Schwarz-ChristofFel map 

Another kind of domains we have been interested in are regular polygons. The conformal 
map applying the interior of the unit circle into the interior of a given polygon is known as 
the Schwarz-Christoffel transformation ^JJI given by 



where Wj are the pre-images of the vertices on the unit circle, ajn are the interior angles at 
the vertices and wq and c are complex constant. A variation of Eq. EZlis the map applying 
the interior of the unit circle on the exterior of a polygon, given by ^Uj 



Finally, the map we are interested in which takes the exterior of the unit circle into the 
exterior of the polygon is obtained by replacing tt; by 1/w in Eq. l38l 

In order to compute these maps, we have used the FORTRAN package SCPACK devel- 
oped by L. N. Trefethen [SZj, as well as the MATLAB SC Toolbox developed by T. A. Driscoll 

m- 

Let us remark that in order to fullfil the condition F{w) w as w ^ oo, we need to 
choose the size of our regular polygons properly. 

4.2.2 Symmetry breaking 

As for the other curves considered, we observe a symmetry breaking of the placement of 
charges on regular polygons. However, the process is slightly different from the one observed 
for the hypotrochoid or the ellipse. Indeed, for both these curves, in the singular limit, the 
symmetry is broken for any finite value of A^. Regular polygons do not present cusps, but 
corners. Depending on the value of the interior angles, the symmetry breaks only for larger 
than a threshold value N* . For small enough values of N, the charges remain far enough from 
the corners to be unaffected by the charges sitting on the other sides of the polygon. When 
N increases, the charges closest to the corner slightly shift to avoid being just in front of the 
charges on the neighboring side. This effect is much more pronounced on curves presenting 
cusps, because the distance between charges placed on two neighboring branches decreases 
much faster as one approaches the cusp. Therefore, the symmetry breaks for any finite N. 

Figure shows the placement of 6 charges on an equilateral triangle and of 8 charges 
on a square. Due to the small value of the interior angles in the case of the triangle, the 
symmetry is broken, while it is preserved for the square. 

4.2.3 Minimum energy configuration 

In order to characterize the minimum energy configuration, we have computed the deviation 
function ^{9) for the regular hexagon, for different values of N. The results are presented in 
Fig. I2UI The charges sitting at the corners remain unshifted, while the ones sitting on one 
side of the hexagon shift toward the closest corner. If the number of charges on one side is 
odd (like for N = 30, N = 36, A = 48 and N = 72), the central charge does not shift. 

For regular polygons, the correlation energy depends on N as with s = 2 according 
to the theory described below. Figures!^ and show this dependence for a regular hexagon 
and a regular octagon. 




(37) 



F{w) = F{wo) + c C"' 11(1 - C/w,r~'dC 



(38) 
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Figure 19: Placement of 6 charges on an equilateral triangle and of 8 charges on a square; 
the symmetry is broken in the first case and preserved in the second. 
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Figure 20: Deviation function -0 for the regular hexagon 
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Figure 21: Dependence of the correlation energy on for the regular hexagon. 
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Figure 22: Dependence of the correlation energy on N for the regular octagon. 
4.3 Scaling arguments for the energy dependence 

The different dependences of the energy on the number of charges N can be easily understood. 

The energy itself is of order unity for a two-dimensional conductor of arbitrary size. There 
are A^^ interactions, each of strength l/N"^ . The size is picked so that this leading order term 
in the correlation energy is zero for when the charges are in their unshifted position. The 
shift in position for each particle is of order and that all the particles and interactions 
participate roughly equally in the shift. Since the energy obeys a variational principle, the 
energy varies as the square of the shift, or as l/N"^. 

To find the shift in the singular case, consider first a simplified model in which charges 
appear in two concentric closely spaced circles. To make all the important distances compa- 
rable take the spacing between the circles to be the same as the spacing between the charges 
in each circle, i.e. 2tt/N. Compare two cases, one with the charges lined up (see Fig. I28|): the 
other with the inner circle rotated so that no two charges are very close to one another (see 
Fig. 124(1 . In the second case, the distance between the closest charges has increased roughly 
by a factor of two. These closest interactions, each of order l/N"^ have changed by an additive 
term of order (1/A^^) In 2. Thus the entire change in correlation energy is of order 1/N . 




Figure 23: Two concentric circles, with the charges lined up. 

The very same argument applies to the hypotrochoid, except in that it applies to only 
some of the charges. Look at the small distance s from the point of the hypotrochoid. There 
are, once gain, two almost parallel lines. At a distance s we have a number of charges, J, 
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Figure 24: Two concentric circles, shifted. 



given by {J/N)"^ ~ s and a separation between the two lines given by (5 ~ (j/N)^. The 
separation between two changes on a single branch is the J derivative of (J/N)'^ or ~ J/N'^. 
The number of charges sensitive to the distortion effects are the ones with 5 less than or of 
the order of this separation. Hence we have 

{J/Nf ~ (J/iV)(l/iV) 

or a number of particle J ~ N^^"^ whose interaction with the nearest neighboring charges can 
be changed by an amount of order unity. These particles each have an interaction of order 
I/N"^, so that the total change is as observed numerically. 

For the very skinny ellipse, or for the line segment, all the particles participate, as happens 
with the concentric circles, so that the change is of order 

For a polygon, only the very few particles nearest the corner participate. Hence the 
singular change in the energy is of order 1/N'^. 

The crossover observed for values of c close to but smaller than c* can be understood as 
follows. For a sufficiently small number of charges, the charges behave as if the curve were 
singular: the symmetry is broken and the above arguments apply. But when N increases, 
the smoothness of the curve prevails, the symmetry is restored and the correlation energy 
goes as A^~* with s = 2. This is illustrated in Fig. I25( where the placement of charges near 
a cusp is shown for the hypotrochoid with c = 0.45, for = 24 and N = 240. 



4.4 Symmetry breaking 

We can also estimate the energy of the symmetric placement of the charges (See Appendix 
IB]) . By considering a test (possibly non-optimal) uniform symmetric configuration Wj = 
exp(27Tij /N+i6o), we see that, with appropriate choices of the minimum energy symmetric 
configuration for an ellipse has an energy 

^ _ logN U[N-^ - N-^]log{l-e'^') N odd, 

" ^ W ~ jiV-a iog(i _ + AT-i log(l - e"^^) even. 

For the hypotrochoid, we have 

^ _ logiV (CN-^e-^' iV/0mod3, 

" ~ W ~ \ C7'iV-2g-iVe/3 + ^,/^-2g-JVe AT = Q mod 3, 
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Figure 25: Charge placement for the hypotrochoid with c = 0.45, for 24 and 240 charges. 
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where we expect that, for an appropriate the constants C,C',C" are ah positive (See 
Appendix |B)). For iV > 1 and e > 0, 

- fog(l - ^ _ i^g(^ _ ^-Ne^^ ^ ^-Ne_ 

Consequently, if N is odd, rather than distributing the charges uniformly on the ellipse, 
a "lower" energy symmetric state would correspond to having [N — l)/2 charges on one 
"branch" of the ellipse (say between = and vr on the unit circle, and having {N + l)/2 
charges on the other branch, with positions corresponding to uniform distribution of — 1 
and N + \ charges respectively. For large iV, the energy of such a "symmetric" state with 
charges at the sharp "ends" of the ellipse will be 

En^-^ + D,N-^ log(l - - D,N-^ log(l - e"^^), 

where Di and D2 > are constants that depend on whether N is even or odd. A similar 
argument for the hypotrochoid yields 

2A 

with a positive constant C that depends on whether or not N is divisible by 3. 
The above arguments show that, the "symmetry broken" states have 

logA^ {an-"" ellipse 

Er, ~ 



2N BN-^I^ hypotrochoid 



with A,B > 0. Comparing the two expressions, we see that the crossover in the energy 
scaling, as well as the symmetry breaking transition occur at 



A-i(l - ellipse, 
'1 



3A-^ (i log N + F) hypotrochoid. 



where F = log B — log C. 

We compare the formula for the hypotrochoid with numerical observations. For the 
hypotrochoid, the above formula yields 



1 

- exp 



31ogA^ 3F 
'2 N ^ aT 



Fig. [2ni is a plot of the c-value for the symmetry breaking, as a function of A, for the 
hypotrochoid. There is one fitting parameter, F, which can depend on whether or not A is 
divisible by 3. We see that the theory agrees very well with the numerics, both for A = mod 3 
and A 7^ mod 3. When A = 0mod3 there are three "corners" which show a distortion from 
the symmetry breaking while in the other cases there is but one. Hence the former case has 
the most to gain from the symmetry breaking, and so we expect the asymmetry to persist 
down to lower values of c for that case. The data supports this expectation. 



4.5 Universality in the Charge locations for singular curves 

We first consider the problem of the placement of A^ identical charges of strength 1 /N on the 
line segment Pi = [—1, 1]. The Fekete polynomial fnix) for Fi is given by 

N 

fN{x) = n 1^ ~ 
3 
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Figure 26: A comparison between the predicted and the numerically observed symmetry 
breaking. 
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where Q are defined by 



TV N 



N N 



nnic-c. 



max 

xi,...,xjve[— 1,1] 



i=l j^i 



i=l jj^i 



/at obeys the equation (HH] 



il-x^)f'j^ + NiN-l)fN = 0. 



The solution can be expressed as 




1 



fN{x) = cn{x^ - 1)PnMx) = 2 



I TV 



{Pn{x) - Pn-2{x)} 



where denotes the Legendre polynomial of degree k j33l 1341 I38j . N identical charges 
therefore place themselves at the ends of the segment and at the positions of the extrema of 
the Legendre polynomial of degree — 1, or equivalently at the intersections of the Legendre 
polynomials of degrees N and N — 2. Figure l271 represents the case = 20. The charges 
accumulate toward the ends. This is consistent with the solution of the continuum problem, 
i.e., the equilibrium charge density on the line segment Fi = [—1,1], corresponding to the 
situation ^ cxd, and given by 



One way to characterize the placement of N identical charges on a line segment is to 
compute the ratio of two successive intervals between the charges, Aj+i/Aj where Aj = 
\xi+i — Xi\, i = 1 corresponding to the charge placed at one extremity of the segment. The 
asymptotic dependence of Aj+i/Aj on N is given by 



where Ri and Ai are two constants depending on i. 

The star symbols in figure 051 show the dependence of Aj+i/Aj on i for the line, obtained 
for N = 200. As expected, Aj4_i/A,; tends to 1 when i increases, i.e., in the central part of 
the line segment. 

The other singular curve we are interested in is the hypotrochoid with c = 0.5. In this 
singular limit, the charge symmetry is broken for any finite number of charges A^. Moreover, 
in the limit of large A^, the placement of the charges near the cusps is similar to the placement 
of charges on the line segment. We have projected the charges positions near one cusp on 
the symmetry axis of this cusp, and we have computed the quantities Aj+i/Aj as defined 
for the line segment. The results are shown in Fig. EH by circles, computed for A^ = 800. 
The correspondence with the results for the line is excellent. This is a signature of universal 
behavior of the charge locations close to the cusp even in the symmetry broken regime. 

Let us remark that the situation is completely different in the case of regular polygons, 
where we consider the placement of charges near a corner instead of a cusp. This is illustrated 




Figure 27: Placement of 20 charges on a line segment. 
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Figure 28: Dependence of Aj+i/Aj on i. 



in Figs. 1291 and l3(Jl which show the placement of 120 charges on a regular hexagon. The second 
figure shows the projection of the charges on the horizontal axis. The shift between two 
facing charges is very small, and their projections remain paired. This behavior is completely 
different from the one observed for the line and the hypotrochoid 




Figure 29: Placement of 120 charges on a regular hexagon. 



OD QD O dD (S) ® O OO* 

Figure 30: Projection on the horizontal axis of the charge positions on two opposite sides of 
the regular hexagon of Fig. 1291 the filled circle is the corner charge. 



5 Discussion 

We have investigated some aspects of the electrostatics of discrete charges in two dimensions. 
The problem is of significant interest because of its connections to approximation theory 
[23 IS] ) constructing conformal maps [30] and Integrable systems |^ . 

A dynamical version of this problem is also closely related to various Laplacian growth 
models because of the deep connections between conformal maps and Laplacian growth 
jl8l 1121 [2|. The dynamics of discrete charges in two dimensions, when they are confined 
by appropriate potentials is therefore relevant to electronic droplet in Quantum-Hall 
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systems 0, and to the Laplacian growth problems including Hele-Shaw and DLA |18j . 
The rich structure of the static problem certainly leads one to expect a similar richness in 
the associated dynamics, and we will investigate this question further. 

In addition, the static problem is also a prototype for the following class of problems - 

Find the "optimal" way to discretize a continuum quantity. 

As stated, this discrete optimization formulation is much too general, for us to say any- 
thing useful about it. Given a specific problem however, there is an appropriate continuum 
quantity, and an appropriate notion of optimality. Such problems arise in many places includ- 
ing approximation theory [HHI , optimal transport jl5| , in discretizing integrable systems 



As in the case of 2-D electrostatics, these discrete problems can have a much richer struc- 
ture than the corresponding continuum problem. In particular, discretizing the continuum 
solution in an obvious manner will not, in general, be a good solution for the discrete problem. 
One explicitly needs to account for the discreteness in the formulation. As we show in Sec. 1, 
the locations of the discrete charges correspond to discretizing a density p which is not the 
continuum solution that minimizes the electrostatic energy. Rather, the density p minimizes 
the modified energy, Es + Ecorr, which are given in Eqs. (|22|) and (|23|) . 

We can increase the "discreteness" of the problem, by making curve more singular, or 
equivalently reducing the total number of charges A'^, for a fixed curve. As we increase 
the "discreteness", we get further away from the continuum solution. In particular, we have 
interesting, qualitatively new phenomena including symmetry breaking, which implies, among 
other things, the lack of a unique minimizing configuration. This is in sharp contrast to the 
continuum problem, which is convex, and therefore has a unique solution. Discreteness can 
thus introduce non-convexity into a problem whose continuum version is convex, and thereby 
change the problem qualitatively EH . 

One of the results from this paper is the universal behavior of the the charge distribution 
, viz., the ip function in the neighborhood of cusps in the symmetric regime. This behavior 
comes out of minimizing the modified energy Eg + Ecorr- This is no longer valid in the 
symmetry broken regime, as illustrated by Figs. [7||H1 and 1^1 Despite this, universal behavior 
persists even in the symmetry broken regime, as we illustrate in Fig. 1281 It is interesting to 
understand the nature of this universality, and we will explore this further in the future. 
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A Saddle point evaluation of Sm 

To evaluate the successive terms Sm in the series solution Eq. (|32j) . we need the following 
result - 

Lemma 1. Let a be a positive real number. Then, 



[HlEni, etc. 




a 



n 



6 \ \m 



) 



) 
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where 



and 



m + n 



(a + V4q! + a'^f 
8a 

(m + n)\ 



^2m — n) {2m — n)\{2n — m)\ 

denotes the binomial coefficient. 

Proof. Setting n = my, replacing the sum over n by an integral with measure mdy, using 
Stirling's formula 



pi = (^^j 



l + O 



we see that the sum 



reduces to 



m 



2m 

E 



[^1 



m 



m + n 



n + mj V 2m — n 



a 



a 



my 



l + y 



(2-y)(2-y)(2y- l)(2y-i) 



dy 



with an O (—) relative error. 

We will evaluate this integral by the saddle point method. Let g{y) be defined as 



9{y) = log 



(2-y){2-y)(2y-l)(2y-i) 



so that the integral reduces to 



m 



V2^ Jy=i/2 ^{l + y){2-y){2y-l) 



dy. 



The maximum for g{y) is at y* given by 

log(a) + log(l + y*) + log(2 - y*) - 21og(2y* - 1) = 0. 
This yields the equation 



ail + y*){2-y*) 



= 1. 



{2y* - 1)2 

Solving this quadratic equation, the relevant solution is 



^13 / a 

2/ = o + o1 



2 2\J A + a 

For < a < oo, we get y* G (1/2, 2), so that the extremum lies within the range of integration. 



lf2y*-lY 1 
exp{g{y )) = - 



a \ 2-y 
Expanding g{y) about y*, we get 

9{y) = 9{y* 



2y/^ 



{a + V^aTc?f 



8a 



a + a — ^/ay 

^ {y-y*f + 0{{y-y*f) 



p. 



2(l + y*)(2-y*)(l-2y* 

so that the point y* is a max;imum for g. We are now in a position to apply the saddle point 
method, and for large m we obtain 



asymptotically, with 0(l/m) corrections. 



□ 
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B The energy of symmetric configurations for singular shapes 

We estimate the energy of the symmetric configuration obtained by distributing the charges 
uniformly on the unit circle, i.e, the points on the curve are zj = F(wj) = i?'[e*(^o+^'^-'/^)]. 
We are especially interested in situations close to the "singular" shape (c = 1 for the ellipse 
and c = 1/2 for the hypotrochoid). For this choice of Wj, we have 



A; / mod N 

J ^ 1 jYe»fc^" A; = mod N 



As we show in Sec. |2 the deviation of the energy from that of a unit circle can be written 

as 

AEr, = ^ [ -^^^InwF'iw). 

J 4:TTIW W" — 1 

Writing this expression in terms of and the matrices a^i and c^, we get 



-ikNdo 



2N 

k>l 1>1 k>l 



We will first evaluate this sum for the ellipse. Using the expressions for aki and Ck in 
Eqs. and respectively, we get 

^^N - 2 2^ —^q,kNe - ^ 2^ —02g,kNe (39) 

g>i '?>i 

We need to do the cases odd and N even separately. 

If is odd, 2q = kN implies that g is a multiple of N and k is even. Therefore 2q = kN 
implies q = jN and k = 2j. Consequently, Eq. yields 

k>i j>i •' 

If N is even, 2g = kN implies that g is a multiple of and there are no restrictions on 
k. Therefore 2q = kN implies q = k(N/2), and A; > 1 is any natural number. Eq (|39|) now 
yields 

fc>l k>l 

We can evaluate these sums using the identity 

°o -Xk 
k=l 

This yields the following exact result for the ellipse 



f \^ - ^1 loefl - e-^(^+2ieo)) Af odd 

^ log(l _ e-^(^+2^^0)) + log(l _ e-iV(e+2ieo)/2) ^ g^g^ ^ ^ 

where, as defined earlier, e = log(c) ~ 1 — c 
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B.l Hypotrochoid 

We can do a similar calculation for the hypotrochoid. Using the expressions in Eqs. H29() and 
((SOI), we obtain 



AE' 



hypo 
N 



1 



Sq+p,kNS2q-p,lN 

p=0 



pj N 



(41) 



We will now look which values of q and p contribute in the above summation. 

q + p = 0,2q — p = mod N =^ 3g = 0, 3p = mod A^. 

We therefore need to do the cases = mod 3 and N ^ mod 3 separately. If 7^ mod 3, 
3q = mod N and 3p = mod imply that that q = jN, p = mN, where j, m are integers. 
In this case, Eq. yields 



1 



^E'r 2Ar 



i>i 



ce 



-31^0 ijAf 



E 

.m=0 



jN 
mN 



1 

iV' 



2^' 



TV 



We can estimate part II of the sum as above. We have 



2Ar2 



E 



ce 



-2^ 



J 



2iV2 



log 



1 - [2ce 



To estimate part /, we use the identity 

N-l 



13=0 m=0 ^ 



Using this, and the above considerations, we get 



[jN 

,mN 



1 

2iV2 



E 



ce 



J 



E0 + 



i>l ■' /3=0 

Interchanging the order of the summations, we get 

N-l 



1 

2iV2 



/3=0 



1 



N' 



Since 



^cos^ 



(l + e^-^/^) =2 
writing Z^^Jo^ = E^=-'^/2' ^i^d using 2c = e"' we get 



7/ 



Af/2-1 



4Ar2 



l3=N/2 



1 + e COS 



^ ■2(-l)/^e-^^cos^(^)cos(3Areo) 



(42) 
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We will estimate this term for 9o = 0. We get, 



^ iV/2-l 

^ E log 



27V2 



l3=N/2 



N 



We will estimate the sums over odd /? and even /? separately. We get 
1 



I = 



2N^ 



a odd 



1 + e-^^cos^ 



7r/3 



+ log 

a even 



l_e-^^cos^ 



TT/? 



Note that the terms with odd (3 are all positive, and the terms with even (3 are all negative. 
It is easy to verify that, for an appropriate eonstant C, we have 

e-''l\l-Ce*)<cos{e)<e-0'l\ 

Using this, for -N/2 < (3 < N/2, we get 



exp 



2N 



1 



iV4 



< COS 



AT 



TTl3 



Tr/2 <e < 7r/2. 



< exp — 



2N 



where C is an appropriate constant. It is clear from the above inequalities, that cos^ 

is significantly different from zero only for < a//V ^ N/2 as TV — oo. Also, if (3 ^ 0{\^) 
we have 



Ar4 



4N Af 



in the region where the exponential e "^^/^^/^^ ig not exponentially small. Consequently, the 
relative error in the approximation 



cos 



N 



7r/5 



exp 



' 2N 



goes to zero as AT — cxd. We thus obtain 



1 

4/V2 



7V/2-1 

E log 

/3=-Af/2 



1 — a cos 



TV 



7r/3 



4iV2 



Ar/2 
-Ar/2 



log 



1 — a exp 



2N 



aVttN 



1 

/ log Tl - Q;exp(-x^/2)l dx, 

kN^J-^ pv / yj , 



where we set x = it(3/\/N, and the relative error is 0{N~^) as iV — oo. 

Note also, if the sum over (3 was restricted to either only odd or only even we have to 
multiply the final integral by 1/2 to obtain the right answer. This can be justified since the 
integrand only varies on a scale (3 ~ Vn 3> 1, so that the sum over the even (3 should very 
nearly equal the sum over the odd /?. 

We will now evaluate the integral in the last line for \a\ < 1. Using 



^ ^ y-^ ^ 



log(l -y) = - 
and integrating the resulting series termwise, we get 

log [l — aexp(— x^/2)] dx 



f 



-Eh 



fc=i 
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We define 



G{a) = ^ 



a" 



k=l 



\a\ < 1. 



It is easy to see that, at a = 1, the sum is the Riemann zeta function G(l) = C(3/2). Also, 
the series for G'{a) diverges like C/\/l — a as a — > 1. Consequently, for a close to 1, we have 

G{a) = CiS/2) - CVT^ + 0{\1 - af^^). 

Also, G(0) = 0,G'(0) = 1 and G"(0) = l/VS. 

Using the above considerations, we can estimate the term I as 



Therefore, we have, 



4\/2iV3 



-Ne 



Ar5/2 



C(3/2) 

g-2JV£ 
'87V372 



+ 



Ne < 1, 
Ne > 1. 



We now consider the case = mod 3. Te only terms in the summations in eq. (|41|) 
that contribute arise from 



q + p = 0,2q — p = mod N 



3g = 0, 3p = mod N. 



Since N = mod 3, 3g = mod N and 3p = mod N imply that that q = jL, p = mL, 
where j, m are integers, and L = N/3. Also, q + p = mod N implies that j + m = mod 3, 
and this also implies 2p — q = pN — {p + q)M = mod N. 

If j = 3k — r, where r € {0, 1,2}, it follows that m = 3t + r, where k, r are integers. Also, 
q = jL = kN — rL and p = mL = tN + rL. Consequently, eq. ()41() yields 



A T^hypo 



Part // gives 



r=0 k>l 
I + 11 



ce 



-mo]kN-rL 



kN -rL 



■fc-[r/2] 

E 

t=0 



kN - rL 
tN + rL 



1 ^ [2ce-3^^o]J'^ 
~2NL^ j 

To estimate part /, we use the identity 

N-l 



2A2 



log 



1 - [2ce-'^*^o] 



2NL ^ ] 



3ieo]N/3 



(43) 



/3=0 



■2ml3r/3 



1 + e 



2ml3/N 



kN-rL 



k-[r/2] 



^ E 



A;A - rL 
tN + rL 



Using this in eq. we get 

2 



1 

2N 



E' 



-27ri/3r/3 



kN -rL 



r=0 k>l 

For all complex number |z| < 1 and r = 0, 1, 2, we have 

■ - log(l - z''), 



E(i+' 

0=0 



(44) 



3E 



z 



3k- 



k>l 



3k — r 



r = 0, 

log(l — z) — uj"^ log(l — — uj log(l — r = 1, 

^ — log(l — z) — log(l — (jjz) — log(l — uj'^z), r = 2, 
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E 

fc>i 



3k— r 



3k — r 



s=0 



where lv = e ^t^/s jg cube root of unity. 

Interchanging the order of the summations in and using the above identities, we get 



N-l 2 2 



13=0 s=0 r=0 



ce ( 1 + e 



2Trif3/N 



N' 



Note that 



Therefore, 



E' 

r=0 



(/3+2s)r 



/3 + 2s 7^ mod 3 
3 s = P mod 3. 



N-l 



2Af2 



/3=0 



1 - e 



-27rj/3/3 



ce 



-34^0 



1 + e 



2-nil3/N 



N 



Breaking this sum up into 6 parts using (3 mod 6, setting = and noticing the similarity 
with the N ^ mod 3 case, we get 



1 

4\/2iV3 



+ 



r=0 



Ar5/2 
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